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We present a theory for a lattice array of weakly coupled one-dimensional ultracold attractive 
Fermi gases (ID 'tubes') with spin imbalance, where strong intratube quantum fluctuations invali- 
date mean field theory. We first construct an effective field theory, which treats spin-charge mixing 
exactly, based on the Bethe ansatz solution of the fD single tube problem. We show that the fD 
Fulde-Ferrel-Larkin-Ovchinnikov (FFLO) state is a two-component Luttinger liquid, and its elemen- 
tary excitations are fractional states carrying both charge and spin. We analyze the instability of 
the fD FFLO state against inter-tube tunneling by renormalization group analysis, and find that 
it flows into either a polarized Fermi liquid or a FFLO superfluid, depending on the magnitude of 
interaction strength and spin imbalance. We obtain the phase diagram of the quasi-fD system and 
further determine the scaling of the superfluid transition temperature with intertube coupling. 



I. INTRODUCTION 

The last few years witnessed an explosion of studies in 
population imbalanced ('polarized') Fermi gases, partly 
fueled by the MIT and Rice experiments with resonantly 
interacting ultracold 6 Li atoms |l|, |2|. One of the motiva- 
tions of these and following experiments [1, 0, S @| is to 
search for the Fulde-Ferrel-Larkin-Ovchinnikov (FFLO) 
state, a spatially modulated superfluid [7fl. So far, no 
evidence of FFLO was found in three dimensional (3D) 
traps. This is in agreement with the mean field theory 
prediction that FFLO state only occurs in a tiny region 
in the phase diagram of 3D polarized Fermi gas @ . 

Contrary to 3D, in one dimension (ID), Bethe ansatz 
@) OH [LH an d Density Matrix Renormalization Group 
(DMRG) studies [H, |l| show that a FFLO-like phase 
(which only has algebraic order) occupies a rather large 
region in the phase diagram of polarized attractive Fermi 
gas. Taking these results together, recently Parish et 
al [3| proposed that the most promising regime to re- 
alize FFLO state is in quasi-one dimension, i.e., a 2D 
lattice array of cold atomic "tubes" [l5| where the trans- 
verse tunneling t±_ is weak. They obtained the phase 
diagram of the quasi-lD system using mean field theory 
[l4T |. These works prompted the experimental search for 
the FFLO phase in quasi-lD atomic gas systems. These 
highly tunable systems are also ideal for investi gat ing the 
lD-to-3D dimensional crossover of superfluids [la ]. 

A few key questions remain open regarding the quan- 
tum phases of quasi-lD atomic systems. 1) What is the 
phase diagram in the limit of small tj_? In this regime, 
intratube quantum fluctuations is dominantly strong [l6| 
so the mean field approach of [l4[ ceases to apply. To find 
the solution, one has to first answer a related question: 2) 
What is the low energy effective description for the ID 
FFLO phase? An effective theory can yield analytical 
expressions for the zero and finite temperature correla- 
tion functions and serve as a starting point to analyze 
the instability of ID FFLO against t± . Note that a field 
theory for ID FFLO state was previously proposed by 
Yang to describe ID superconductors in magnetic field 
[TtJ ■ In Yang's model, spin and charge degrees of free- 



dom are assumed decoupled. However, we shall show 
that this model does not apply to attractive Fermi gases 
with finite spin imbalance ('magnetization') because the 
assumption of spin-charge separation breaks down [l3l ]. 

In this paper, we answer these open questions using a 
combination of several well established theoretical tech- 
niques. We first develop an effective field theory for the 
ID FFLO state, which properly treats the spin-charge 
mixing in trapped two component atomic gases and fully 
takes into account the intratube quantum fluctuations. 
Then we treat the transverse tunneling perturbatively 
[Tq] to obtain the phase diagram of a 2D lattice array of 
weakly coupled ID Fermi gases, which is an anisotropic 
3D system. We also derive the scaling formula for the 
transition temperature of the quasi-lD FFLO state. 

Our approach, based on exact solutions, is complemen- 
tary to the mean field method of Ref 14], which applies 
to the opposite limit of large t±. Taken together, they 
can provide a more complete understanding of the di- 
mensional crossover of quasi- ID attractive Fermi gases 
with finite population imbalance, and provide guidance 
for the experimental search of the quasi-lD FFLO phase. 

II. MODEL 

We consider a gas of two-component fermionic atoms, 
e.g., 6 Li atoms in hyperfine states \F ~ 1/2, to/ — ±1/2) 
(labeled with spin up and down respectively), confined in 
a highly elongated one dimensional (ID) trap interacting 
via an attractive contact potential, 

h 2 N d 2 Nr Nl 

i— 1 1 i— 1 j — 1 

This Hamiltonian is known as the Gaudin-Yang model 
[11, El- Here m is the fermion mass, and the attractive 
interaction stength g = —2h 2 /maiD < 0, a\D is the ef- 
fective ID scattering length and related to the 3D s-wave 
scattering length a s by [2(| HH 

a 1D = (1.03 - — )a±, 
a s 
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where a± = yj fi/mwi is the oscillator length of the trans- 
verse confinement of frequency u>± w x . Let N-\ and 
Ni be the number of spin up and down particles respec- 
tively, and N = JVf + N± the total particle number. We 
define the linear particle density = N/L, the popula- 
tion imbalance p — (JVj — Ni)/N, and the dimensionless 
interaction parameter 7 = 2/ {jiliid)- These are experi- 
mentally controllable parameters. Some of these numbers 
are listed in Ref. [22|. In the above model Hamiltonian, 
we have neglected the slowly varying Gaussian trap po- 
tential in x direction, V(x) — mu> 2 x 2 /2, and assumed 
periodic boundary condition for the system of length L, 
as we are interested in the collective behaviors of the pop- 
ulation imbalanced gas in the thermodynamic limit. The 
effect of the trap can be subsequently taken into account 
using local density approximation. 

For technical convenience, we also consider a closely 
related lattice model, the ID attractive Hubbard model, 



H v = -^(ct CT Ci +1>CT + h.c.) ' 



(2) 



with linear particle density nj, = n/a and population im- 
balance p. Here, n — rij + n± is the particle number per 
site (band filling), and a is the lattice spacing. We define 
the Fermi wave vector fc/ CT = nn^/a for each spin species 
(7 =|, I, and the FFLO wave vector q+ = fe/f — kfi- The 
Hubbard model can be viewed as the (regularized) lat- 
tice version of the continuum Gaudin-Yang model, with 
a playing the role of the short distance cutoff. Explicitly, 
the Gaudin-Yang model corresponds to Hjj in the contin- 
uum limit, a — > 0, for fixed particle numbers and system 
length L, i.e., n <C 1. In this limit, the parameters of 
Hqy are related to those of H\j by [23| 
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Finally, the experimental setup [24j | searching for quasi- 
ID FFLO state consists of a 2D square lattice array of 
ID tubes, each described by Hqy and coupled to its 
nearest neighbors by transverse tunneling of amplitude 

t±, 



H± = -t 



dx V [4, v Jx)i>i,vA x ) + h - c -\- ( 4 ) 



Here ipu=B,/L,cr is the fermion field operator for right/left 
movers with spin a =|,I (25[, and labels nearest 

neighbor tubes. Also relevant to our discussion is the 
intertube Josephson coupling of the form 



Hj = -J I dxJ2i A H x ) A j(x)+h.c.]. 

(id) 



(5) 



The pair operator is defined as A(x) = tj)R,\ (x)ipz,l(x), 
and the pair tunneling, J oc tj_, is generated by sec- 
ond order hopping process. Intertube coupling estab- 
lishes transverse coherence between tubes and drives the 
system through a lD-to-3D dimensional crossover [l6| . 



III. SPIN-CHARGE MIXING IN AN 
IMBALANCED FERMI GAS 

We first obtain a low energy effective theory for the 
Gaudin-Yang model using standard Abelian bosonization 
[25l [26j . which is valid for weak interaction (small 7). 
We linearize the spectrum of each spin species around its 
Fermi points ±fcj CT . Note that the two Fermi momenta 
differ by 6k f — — kfi, since 7^ N±. We also 
define the average and difference of Fermi velocities, «/ = 
(vfi+Vfi)/2 and 6vf = U/f — v/x- We follow the notation 
convention of Ref. [2 a ] , where the bosonization identity 
takes the form 

M x) = ^e ik f^e-^-^, 
yZira 

U L . a 



Here 9 a is the dual field of boson field </>„, Ur/l^ are 
the Klein factors, and a is the short distance cutoff. The 
bosonized Gaudin-Yang model |T]) becomes, 



dx 
2^ 



v f {Ve c ) 2 + v+(y<j) c f + 6v f 6k f 



+ v S (ye s ) 2 + vj(V4> s ) 2 + V2vj6k f vct> s 



V0c 

V2 



6vf(V6 c V9 s + V0 c V0 s 



i(V80 s ) . (7) 



Here the charge and spin field <p c / s = (0j ± 4>{)/y/2 are 
defined in the standard way. We have introduced shifted 
spin fields 

4>s = 4>s - 6kfx/V2, S = 9 S , 

and a short hand notation = Vf ± g/ir. 

For the unpolarized case, p = 0, both 6k and 6vf van- 
ish, so Hb is reduced to the sine-Gordon model. With a 
negative coefficient (g < 0), the cosine term is relevant. 
It is well known that the system, an attractive Fermi gas, 
is a Luther-Emery liquid where spin and charge separa- 
tion holds: the charge sector is gapless and described by 
a Gaussian Hamiltonian and the spin sector is gapped. 
The algebraic decay of the s-wave paring susceptibility, 
which is most diverging, is determined by the Luttinger 
parameter in the charge sector. 

Finite population imbalance brings several significant 
differences. First, it introduces a linear "source" term for 
the (shifted) spin field which acts as an effective magnetic 
field, 

h = V2vj6kf. 

But more importantly, it leads to coupling between the 
spin and charge sector, described by the density and cur- 
rent interaction terms (y9 c V9 s and V(/) c V0 s ) of order 
6vf in Hb- Therefore, Eq. ([7]) clearly shows that a ID 
attractive Fermi gas with finite p is in general not a spin- 
charge separated Luther-Emery liquid. 
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In the limit of very small p, both Svf and 5k / are small 
so we can neglect high order terms oc 5vf5kf in Eq. ([7]). 
Furthermore, in situations where we are allowed to ne- 
glect the spin-charge mixing terms oc Svf but keep the 
linear V<^ s term oc 5k f, then Eq. ([7]) reduces to the spin- 
charge separated model proposed by Yang to describe 
quasi-lD superconductors in magnetic field [l7j ]. 



v f {ve c ) 2 + v+{V4> c f 



+ v f (Ve s ) 2 + vJ(V4> s ) 2 + + cos 



Its charge sector is a massless Gaussian, while its spin 
sector (the second line) is a sine- Gordon Hamiltonian in 
an effective magnetic field h, which we call H s ch- H a Gh 
has been investigated thoroughly by many authors since 
the pioneer works of Ref. [27|, [28|. Increasing h over a 
critical value h c will trigger a phase transition from a 
Luther-Emery liquid (BCS-like superfluid with algebraic 
order and spin gap) to a FFLO-like phase in which the 
spin excitation is also massless and the pair correlation 
function oscillates at wave vector q+ on top of the alge- 
braic decay. The spin sector of the ID FFLO phase of 
Hy can be described as a Luttinger liquid of spin soli- 
tons [13] • These solitons are the ID counterpart of the 
static order parameter domain walls (more precisely, the 
Andreev bound states at the domain walls populated by 
the majority spins [141 ]) which form a crystal in high di- 
mensional FFLO phases [29| , 

Based the spin-charge separated Hamiltonian Hy, 
Yang also pointed out that the magnetic field driven 
BCS-FFLO transition is in the universality class of 
2-dimensional classical commensurate-incommensurate 
transitions [13] ■ To see this clearly, we can solve H s ch 
analytically along the Luther-Emery line to map it onto 
a fr ee ( massive) fermion Hamiltonian by refermionization 
[IS>IIg| : the magnetic field h then plays the role of chemi- 
cal potential, while g controls the size of the band gap. As 
h is increased to h c , the Zeeman energy exceeds the band 
gap, so the upper band becomes populated by spinless 
fermions (spin solitons) . This general picture holds even 
away from the Luther-Emery line, where spin solitons 
become interacting with each other. Thus, Yang's model 
Hy gives an appealing physical picture for the ID FFLO 
like phase in terms of solition liquids, provided that the 
basic assumption of spin-charge separation holds. 

Yet, the spin-charge mixing terms in Hg have scal- 
ing dimensions of 2 and are marginal operators in the 
renormalization group sense. So they cannot be simply 
neglected even if their values are small. Moreover, in cold 
atom experiments the population imbalance < p < 1 
can be arbitrary and usually not small. This is drasti- 
cally different from typical situations in solids where the 
Zeeman energy is much smaller than the Fermi energy. In 
general, the spin-charge mixing terms have to be treated 
on the same footing as other terms in Hb- 

While a direct perturbative analysis of Hb is possible, 
it is unclear how well the bosonized Hamiltonian Hb de- 



scribes the system for the case of strong interaction and 
large population imbalance. For these reasons, we shall 
adopt an alternative but much more powerful method, 
which is valid for arbitrary (n,p,g) and treats the spin- 
charge mixing effect exactly, to construct the bosonized 
effective field theory of Hoy- The method is based on 
the Bethe ansatz solution of the microscopic model Hgy 
(or Hjj) and general principles of conformal field theory. 



IV. EFFECTIVE THEORY FOR THE ID FFLO 
PHASE 

In this section, we show that despite the spin-charge 
mixing, the ID FFLO phase features a more general de- 
coupling of two critical degrees of freedom. In the lan- 
guage of conformal field theory, the ID FFLO state is a 
critical phase with two gapless normal modes and it is 
formally described by the direct product of two Virasoro 
algebras with central charge c — 1 [30]. Such a "two 
fluid" description is in accordance to our physical intu- 
ition in higher dimensions, e.g., the 3D attractive Fermi 
gas with population imbalance in the strong interaction 
(BEC) limit can be described by a fermion-boson mixture 
with some residue interactions. However, the distinction 
between bosons and fermions is lost in ID. Our task is to 
identify the two normal modes (elementary excitations) 
and establish their relation with the usual spin/charge 
excitations. 



A. Summary of main results 

We first present the final results, so readers not inter- 
ested in technical details may skip the derivation and pro- 
ceed directly to the next section. The ID FFLO state is 
a two component Luttinger liquid described in bosonized 
form by the effective Hamiltonian 



-ffFFLO — ^ 

i=l,2 



dx 
2^ 



ui [(vtfO 2 + (v^) 2 



(8) 



The two normal modes (i = 1, 2) are decoupled and have 
different group velocities, u\ and ui- Each mode, de- 
scribed by the boson field ipt and its dual is a super- 
position of the spin up and down fields <j) a and 9 a , 



Ts-l 



<f>i 



0i 

#2 



(9) 

The superscript T means matrix transpose, and the so- 
called dressed charge matrix is given by 



Z = 



Z C c 



Z S c 



z s . 
-z R 



(10) 



Eq. (JH|) clearly shows that in general each normal mode 
ipi is a mixture of the spin field <f> s and the charge field <f> c . 
Therefore, elementary excitations of the FFLO state are 
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fractional states carrying both charge and spin, i.e., the 
superposition of spinons and holons [3lJ . The degree of 
spin-charge mixing is determined by the dressed charge 
matrix Z, or equivalently matrix Y = [Z T ] _1 (30l.[3ll[33. 
l33l ] . As shown in Ref. [3l[ , spin-charge separation occurs 
only if the charge matrix 



.9 



Z _1 Y 



possesses a Z 2 symmetry, g n = g 22 , 512 = 521, i.e., g is 
invariant under the exchange 1 <-» 2. All the parameters 
in the effective bosonized Hamiltonian Hfflo, such as 
{Zij,Ui}, can be computed from the parameters of the 
microscopic Hamiltonian Hgy or Hjj. Eq. ([8|)- (fT0)l are 
the main results of this section, they form the basis for 
all subsequent discussions on correlation functions and 
quasi- ID phases. 



B. The method 

Below we outline the main steps leading to the effec- 
tive field theory Eq. (jHJ) and the algorithm to compute 
parameters {Zij,Ui\. The procedure takes a detour to 
make use of several well established, nontrivial results 
for the ID Hubbard model. 

(I) We start from the attractive Hubbard model Hjj- 
Its continuum limit, summarized in Eq. ([3]), yields the 
Gaudin-Yang model Hgy- Therefore, in this limit, the 
effective field theory of the two models are the same. 
Indeed, there is a one to one correspondence between the 
Bethe ansatz integral equation for the Hubbard model 
[23l | and the Gaudin integral equation for Gaudin-Yang 
model 34]. In principal, all analysis can be performed 
on the Gaudin-Yang model directly without resorting to 
its lattice regularization to arrive at the same conclusion. 
In practice, however, the latter route is less economical 
since the Bethe ansatz solution of the continuum model, 
especially its dressed charge formalism, is less developed. 

(II) We map the attractive Hubbard model Hjj for- 
given (U/t,n,p) to a repulsive Hubbard model H' v with 
onsite repulsion U, density n' = 1 — np, and imbalance 
p' = (1— n)/(l— np) bythe well-known staggered particle- 
hole transformation [23, K3a , 



"id- 



(11) 



The motivation behind this is two- fold. First, it is tech- 
nically more convenient to solve H'jj via Bethe ansatz 
f35l ]. Second, as shown by Penc and Solyom [38| based 
on the work of Frahm and Korepin [301 ] . the repulsive 
Hubbard model at general (non-half) filling and mag- 
netic field is equivalent to a two-component Tomonaga- 
Luttinger (TL) liquid. We shall exploit this result to 
construct the field theory for Hjj (and Hgy) from the 
Bethe ansatz solution of H' v . 

(Ill) For given (U/t,n',p'), we numerically solve the 
coupled Bethe ansatz integral equations for H[j. We 
first determine the density distribution functions p(k) 



and c(A), defined within the Fermi interval k £ (—Q,Q) 
and spin rapidity range A 6 (—A, A) respectively. This 
requires solving 

1 f A 

p(k) = hcosfc / dAai(sinfc — A)er(A), 

27r J -A 

rQ 

cr(A) = / dkp{k)a\(sva.k — A) 
J-Q 

- f d\'a(\')a 2 {\-\'), 

J -A 

where a n (x) = nlJ / {Airi) / [{U / At) 2 + x 2 ], together with 
the number constraint equations, 

rQ 

/ dkp(k) = n , 
J-Q 

J A d\a{X)=n\ = ^{l-p'). 



This is done by an iteration procedure which also fixes 
Q and A self-consistently. Each integral equation is con- 
verted to an algebraic equation using Gaussian-Legendre 
quadrature. The code was checked to yield exactly the 
same result as those in Ref. [13, HEJ . 

Once Q and A are found, it is then straightforward to 
find other quantities: the dressed energies [e c (/c), e s (A)], 
their derivatives, [e' c (k), e' s (X)], and the dressed charges, 
[£ cc (fc),£ cs (A),£ sc (fc),£ss(A)], by solving corresponding 
integral equations. These lengthy equations are well doc- 
umented and will not be repeated here [H,[30, 38]. The 
dressed charge matrix elements are given by 

Z C c = ^cc(Q), Z cs = ^ CS (A), 

z sc = u(Q), z ss = £ ss (A). 

The normal mode group velocities are, 



"1 



e' c (Q) 
2irp(Q) 



U 2 



<{A) 
2tt<j(A) ' 



Note that the subscripts c and s here and in the Bethe 
ansatz literature have to be distinguished from the charge 
and spin index in Abelian bosonization such as in Eq. (JT)). 
This fact has been emphasized for example by Voit et al 
[39l ] . To avoid confusion, we adopt index i = 1,2 to label 
the two normal modes instead of c and s. 

(IV) We map H^ to a two-component Tomonaga- 
Luttinger Hamiltonian, H' TL , following Penc and Solyom 
(3^ | . It becomes transparent in bosonized form 



■<74,T -52,T)( V6, t) 5 



+ ( v fA + .94,T + .92,t)(V0 T ) 2 
+(vf,l + 54,1 -52,l)(V0i) 2 
+ 54,1 +52,l)(V0i) 2 

2(ff 4 ,j_ - 32,jJV0 t V9 l + 2( 54 ,± + .g 2 ,j_)V0 T V0x] . 
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The Luttinger parameters in the g-ology notation, g, are 
given in terms of {(/,-. Z, j) by Eq. (G.12) to (6.17) in Ref. 
l38l . For example, 

Vf,-[ + ff4,T - 92, T = u i( z cc - Z sc ) 2 + u 2 (Z cs - Z ss ) 2 , 

, - , - Ul^+tta^ 

v fA + 94,1 + 92,1 = T^-^ 7 7 \2 - 

{Access — ^cs^sc) 

The low energy effective theory H' TL yields the same fi- 
nite size spectrum and correlation functions as the orig- 
inal microscopic model H' v . It is instructive to compare 
this nonperturbative result with that obtained by weak 
coupling bosonization, Hb in Eq., 0. 

(V) Finally we particle-hole transform H' TL back to ob- 
tain the effective Hamiltonian Htl for Hjj. In bosonizcd 
form, Htl is obtained from H' TL by replacement <f>^ — > 
— fa and #j — > — #4. Since Htl is quadratic, we diag- 
onalize it by a canonical transformation using matrix Z 
[Ulill. The end result is Eq. ©. 

V. SINGLE PARTICLE AND PAIR 
CORRELATION FUNCTIONS 

The bosonized effective Hamiltonian Eq. ((5J) makes it 
straightforward to compute the correlation function of 
any operator of the form 

0(x) = e^k^+^kx)] 

where a and b are arbitrary two-component vectors, and 
</> T = 0J,), § T = (9 h 0i)- Using the relation Eq. ©, 
we find that the correlation function (O(x)O'(0)) has 
the asymptotic behavior of x~ 2S at zero temperature as 
x — > oo. The scaling dimension of O is given by 

6=~[a T g- 1 a + b T gb}. 

The charge matrix is defined by g = Z Y. 

We apply this result to compute the single particle 
propagator and the singlet pair susceptibility 

G a {x,r) = -{T T i> Rt(T (x,T)^ Ri(T (0,0)), 
Xi x,r) =-(T T A(x,r)At(0,0)). (12) 

We shall focus on G4 for the majority spin (f), since G4 
only has subdominant power law divergence (its scaling 
dimension is larger than that of G4, i.e., 5^ > <5j). For 
example, from the bosonized representation of fermion 
field operator i/jr^, Eq. ([6]), we find a — (— 1,0) T and 
b = (1,0) T for G4. It follows that for x — * 00 at zero 
temperature, 

G r(*)~^p <5T=ft+&, (13) 

Ml' <^A = 7l+72- (14) 



The scaling dimensions of the single particle field oper- 
ator and pair operator A are directly related to the 
dress charge matrix elements that describe the effect of 
spin-charge mixing, 

Pi = \[{Zuf + (Yii) 2 }, 02 = \[{Z2i? + (Y 21 n 

1 Z 2 1 Z 2 

71 = l [Zl + Jd^W h 72 = l [Zl + fdet% ] - 

(15) 

Here Z\\ denotes the (1,1) element of the matrix Z. We 
also define the corresponding anomalous dimensions 774 — 
2 - 2S h and tj a = 2 - 25 A . 

Eq. (fMj) is the hallmark of ID FFLO phase: the s-wave 
pair correlation function oscillates in space at the FFLO 
wave vector q+ = fc/j — on top of the algebraic decay. 
Although it looks similar to Yang's conclusion based the 
spin-charge separated model Hy , the present result fully 
takes into account the effect of spin-charge mixing and 
remains valid for all values of < p < 1 . Our results can 
also shed more light on the nature of the magnetic field 
driven BCS-FFLO transition in ID. Under a particle-hole 
transformation, it is equivalent to the Mott transition of 
repulsive Fermi gas in magnetic field (driven by varying 
chemical potential). As shown recently by Frahm and 
Vekua [32j, due to the intrinsic spin-charge coupling, it 
is not in the same universality class as the single mode 
commensurate-incommensuarate transition. 

At finite temperature T, G|(x,r) and xi x i T ) still 
factorize due to the decoupling of two normal modes. 
After some algebra, we obtain the Fourier component 
Go = G|(/c = fc/.|, u> = 0) (up to a phase factor), 

Go = (2^ 3 aT 2 )- 1 I 1 (^) TT C—Y^V^- (16) 

112 u * 

Similarly, xo = x(^ = 9*,^ = 0) is found to be 

Xo = (2n 2 aT)- 2 l 2 (^) J] (^) 2 ^. (17) 

Here, the dimensionless function Ii/ 2 are denned as 

/OO r7T 
dx / ^sinirj[sinh 2 (xy i_ i) + sin 2 (t)]" ft , 
-00 Jo i 

/OO />7T 
dx / dtTTlsinh^^-^ + sin 2 ^)]- 7 *. 
-00 Jo i 

These results will become useful in the RPA calculation 
for the quasi-lD systems. Note that the ratio u\jui is 
not a rapidly changing function of 7. Fig. [T] shows its 
value for p = 0.2. 

VI. PHASE DIAGRAM OF THE QUASI-1D 

SYSTEM AT T = 0. 

We now apply the effective theory for single tube to 
study the quasi-lD system realized in experiments: a 
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FIG. 1: The ratio U1/U2 as function of 7 for p = 0.2. 

square lattice array of such tubes coupled by transverse 
hopping t±. We are interested in the fate of the ID FFLO 
phase as the t± is turned on. 

In the limit of small t±, we can treat H± and Hj as 
perturbations to the individual ID Hamiltonian -//fflo- 
In the language of renormalization group (RG) , both H± 
and Hj are relevant perturbations: the scaling dimen- 
sions of the single particle and pair tunneling arc both 
smaller than 2. The first-order RG equations for the ef- 
fective inter-tube couplings tx(K) and J(k) at momentum 
scale k (with k — > 0) read: 

Kdtx{K)/dK = (2Sf - 2)tx(/c), ( .o\ 
KdJ(n)/dK = (25 A - 2)J(k), 1 j 

where the equation of fj_ is for the spin f tunneling (as 
we have emphasized before, it is more relevant than the 
spin j). In other words, the fate of the ID FFLO phase 
is controlled by the relative magnitude <5j and 6 a [l7l ]. 
For <5a < <5fi pair tunneling is most relevant and the 
system flows into a quasi-lD FFLO state. In this state, 
strong effective Josephson coupling locks the phases of all 
tubes to establish the overall phase coherence to produce 
a genuine superfluid state. For 5 a > <5|, however, single 
particle tunneling is most relevant and the system flows 
into a partially polarized Fermi liquid (FL) state with 
well defined quasiparticles. In the latter case, the actual 
ground state of the quasi- ID system at zero temperature 
depends on the residue interactions between quasiparti- 
cles, the details of which are not captured by the leading 
order RG analysis presented here [17| . Such limited pre- 
dictive power is inherent to all leading order RG analysis 
on coupled Luttinger liquids. Therefore, the Fermi liquid 
state predicted here should be understood as a region in 
the phase diagram where the superfluid transition tem- 
perature is significantly suppressed by the weakening of 
effective intertube Josephson coupling. It is important 
to bear in mind that other instabilities may take over at 
lower temperatures leading to a ground state with broken 
symmetry. 




FIG. 2: Zero temperature phase diagram of a quasi- ID at- 
tractive Fermi gas in the limit of weak inter-tube tunneling. 
7 is the interaction strength of the Gaudin-Yang model and 
p is the population imbalance. 



The T = phase diagram of the quasi- ID gas system 
based on leading order RG is shown in Fig. [2] It is ob- 
tained by going through the steps outlined in Sec lIVBI 
and taking the continuum limit, a — > 0, for fixed parti- 
cle numbers and system length L [see Eq. As the 
intertube tunneling is turned on, the ID FFLO phase 
(originally occupying the whole region < p < 1 for all 
7) splits into two distinct phases, an FFLO superfluid 
and a polarized Fermi liquid (FL). Intuitively, stronger 
attractive interaction (larger 7) favors the FFLO phase. 
From Fig. [5J one can read off the critical interaction 
strength required to realize the quasi-lD FFLO state for 
given p. For fixed interaction 7, increasing imbalance 
would drive the system out of FFLO into a Fermi liquid 
phase. A crucial feature of the phase diagram is that the 
FFLO phase survives in a smaller region in quasi- ID than 
true ID (single tube). This shrinking trend observed at 
ix — ^ is expected to continue as tj_ increases, since we 
know that FFLO state in 3D only occupies a tiny part 
of the phase diagram y|. Fig. [3] shows T = phase di- 
agram of the quasi-lD attractive Hubbard model in the 
limit of t j_ — > 0. We observe a similar splitting of the 
ID FFLO phase. Ref. |37| also discussed the phases of 
weakly coupled Hubbard chains in the presence of finite 
spin polarization, close in spirit with ours. While our 
main focus here is the continuous gas systems (without 
lattice in the x direction), our result of quasi- ID Hubbard 
model (Figure 3) agrees with Ref. |37|. 



VII. THE TRANSITION TEMPERATURE OF 
QUASI-1D FFLO STATE 

In this section, we use the random phase approxima- 
tion (RPA) [H, EH to compute the FFLO superfluid 
transition temperature T c for the quasi-ID system with 
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FIG. 3: (Color online) Zero temperature phase diagram of 
the quasi- ID attractive Hubbard model, Hu + Hj_ + Hj, in 
the limit of small t±. p is the population imbalance. Each 
line for fixed density n (top panel) or interaction strength U/t 
(bottom panel) is the phase boundary separating the FFLO 
and Fermi liquid (FL) phase (indicated by arrows). Inset: 
phase diagram of the ID attractive Hubbard model (single 
tube). 



a small but finite t±, going beyond the perturbative RG 
analysis in the previous section. RPA has been success- 
fully applied to study coupled Luttinger liquids (p = 0) 
[4lj . Here we generalize it to the case of p > with 
spin-charge mixing. Within RPA, the 3D single particle 
propagator 

9f 1 (k,k± = 0,w n ) = G^(k,u n ) - z±t±, (19) 
and the 3D pair susceptibility 

t r- 1 (^,fcx=0,cj = 0)=x c 7 1 -^^ (20) 

where k± is the transverse momentum and z± = 4 is the 
transverse coordination number of every single atomic 
gas tube in a square lattice array. The ID propagator 
Gj, and susceptibility xo have been computed in Section 
|V]and are given by Eq. (fT!)|) and (fTT)) . respectively. 




FIG. 4: RPA predictions for the t± dependence of the single 
particle crossover (Luttinger liquid to Fermi liquid) temper- 
ature, Tfl ~ ij_ ^ , and the FFLO transition temperature, 
T c ~ t 2 / 71A , as functions of 7 for p = 0.2. 



RPA predicts two temperature scales that character- 
ize the lD-to-3D crossover as the temperature is lowered. 
^-[■(fc/T, 0, oj n =i) starts to develop a pole at the single par- 
ticle crossover temperature Tfl- From Eq. p9[) . we find 



- FL 



w 



tj_Z A 



W 



Ii( 



Ui , 



1 /V7 



n 



(l/2-2ft)/»?T 



(21) 

Similarly, the divergence of the 3D pair susceptibility 
^T(q*,0, 0) defines the two particle crossover tempera- 
ture T r , 



Jzj_ Ul 

L 2{ ) 

Vf Zi 2 



1 l/'JA 



(l/2-2 7i )/7,A 



(22) 



Here the ultraviolet cutoff W — Vf/a is roughly the band 
width (or Fermi energy). Above max^FL, T c ), ther- 
mal fluctuation destroys coherence between the tubes, 
the system behaves as uncoupled two-component Lut- 
tinger liquids with fractional excitations, each described 
by Hfflo- As the temperature is lowered, if Tfl > T c , 
these fractional excitations first confine into sharply de- 
fined quasiparticles at temperature scale Tfl, and the 
system crosses over to a partially polarized Fermi liquid 
(the Fermi liquid may develop other instabilities at lower 
temperature depending on the details of the residue inter- 
actions between the quasiparticles). On the other hand, 
if T c > Tfl , the system first undergoes a phase transition 
into an FFLO phase with long range order. In this case, 
T c can be identified as the superfluid transition temper- 
ature. In the limit of vanishing p, the transition into the 
FFLO superfluid can be viewed as the condensation of 
soliton liquid into crystal, i.e. static domain walls (l7l |. 

Whether the single or two particle crossover occurs first 
depends on the interaction strength. At low densities, 
we find from Eq. (ST) and (22) that Tfl > T c at small 
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7 (weak interaction) while T c > Tfl for large 7 (strong 
interaction). Thus, the schematic plot of Tfl and T c as 
functions of 7 looks similar to Fig. 1 of Ref. [H[ for 
coupled Luttinger liquids. At some critical value of 7, 
Tfl and T c are comparable to each other. This marks the 
phase boundary between the Fermi liquid and the FFLO 
phase. All these results are in qualitative agreement with 
those obtained in the previous section from RG analysis. 

Eq. (|22| shows that T c scales with J as a power law, 
with the exponent given by the inverse of anomalous di- 
mension r)A, 

T c cx J 1 /"*. 

The exponents 77^ 1 and r)7 are plotted in Fig. [4] for p = 
0.2. We observe that at weak interaction Tfl oc t±, while 
at strong interaction T c cx J cx tj_ . The growth of T c with 



inter-tube coupling will eventually stop when t± becomes 
so large that it can no longer be treated as a perturbation. 
Then the system becomes more 3D-like, and T c starts to 
drop as t± is further increased [14j . Our results support 
the argument of Ref. [LjJ that the optimal value of T c is 
realized for small but finite t±. 
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